!-------------------------------------------------------------LICENSE--------------------------------------------------------------!
!                                                                                                                                  !
!The MAP code is written in Fortran language for magnetohydrodynamics (MHD) calculation with the adaptive mesh refinement (AMR)    !
!and Message Passing Interface (MPI) parallelization.                                                                              !
!                                                                                                                                  !
!Copyright (C) 2012                                                                                                                !
!Ronglin Jiang                                                                                                                     !
!rljiang@ssc.net.cn                                                                                                                !
!585 Guoshoujing Road. Pudong, Shanghai, P.R.C. 201203                                                                             !
!                                                                                                                                  !
!This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License         !
!as published by the Free Software Foundation; either version 2 of the License, or (at your option) any later version.             !
!                                                                                                                                  !
!This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of    !
!MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more details.                        !
!                                                                                                                                  !
!You should have received a copy of the GNU General Public License along with this program; if not, write to the Free Software     !
!Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.                                                   !
!                                                                                                                                  !
!-------------------------------------------------------------LICENSE--------------------------------------------------------------!

!==================================================================================================================================|
subroutine write_silo(file_unit, var_name, var, x, y, z, nx, ny, nz, n_var, var_dim, block_id_char)
!==================================================================================================================================|

   implicit none
   include 'silo.inc'

   integer(4), parameter :: n_dim = 3
   integer(4), parameter :: n_point_cell = 8

   character(*), intent(in) :: var_name(n_var)
   character(256), intent(in) :: block_id_char
   integer(4), intent(in) :: nx, ny, nz, n_var, file_unit
   integer(4), intent(in) :: var_dim(n_var)
   real(4), intent(in) :: x(nx), y(ny), z(nz)
   real(4), intent(in) :: var(sum(var_dim), nx, ny, nz)

   character(256) :: mesh_name_block
   integer(4) :: error, silo_status, index_var, i_var

!----------------------------------------------------------------------------------------------------------------------------------|
!  prepare
!----------------------------------------------------------------------------------------------------------------------------------|
   error = dbmkdir(file_unit, trim(block_id_char), len_trim(block_id_char), silo_status)
   error = dbsetdir(file_unit, trim(block_id_char), len_trim(block_id_char))

   mesh_name_block = 'amr_mesh'
   error = dbputqm(file_unit, trim(mesh_name_block), len_trim(mesh_name_block),                                                    &
               'x', 1, 'y', 1, 'z', 1, x, y, z, (/ nx, ny, nz /), n_dim, db_float, db_collinear, db_f77null, silo_status)

   error = dbset2dstrlen(8)
   index_var = 1
   do i_var = 1, n_var
!print*, trim(var_name(i_var)), len_trim(var_name(i_var))
!----------------------------------------------------------------------------------------------------------------------------------|
!     write the scalar data
!----------------------------------------------------------------------------------------------------------------------------------|
      if (var_dim(i_var) .eq. 1) then
         error = dbputqv1(file_unit, trim(var_name(i_var)), len_trim(var_name(i_var)),                                             &
            trim(mesh_name_block), len_trim(mesh_name_block), var(index_var, :, :, :), (/ nx, ny, nz /), n_dim,                    &
            db_f77null, 0, db_float, db_nodecent, db_f77null, silo_status)
         index_var = index_var + var_dim(i_var)
!----------------------------------------------------------------------------------------------------------------------------------|
!     write the vector data
!----------------------------------------------------------------------------------------------------------------------------------|
      elseif (var_dim(i_var) .eq. 3) then
         error = dbputqv(file_unit, trim(var_name(i_var)), len_trim(var_name(i_var)),                                              &
            trim(mesh_name_block), len_trim(mesh_name_block), n_dim,                                                               &
            (/ "vector_x", "vector_y", "vector_z" /), (/ 8, 8, 8 /), var(index_var : index_var + 2, :, :, :),                      &
            (/ nx, ny, nz /), n_dim, db_f77null, 0, db_float, db_nodecent, db_f77null, silo_status)
         index_var = index_var + var_dim(i_var)
      endif
   enddo

   error = dbsetdir(file_unit, '..', 2)

!----------------------------------------------------------------------------------------------------------------------------------|

   return
end subroutine write_silo
